#include <stdlib.h>
#include <stdio.h>
#include "slave.h"
#include "3d7p_aceMesh.h"
#include "task.h"
#include "ldm_malloc.h"
 __thread_local volatile unsigned long get_reply,put_reply;
//kernel region
__thread_local_fix double *A_lhs_array = NULL;
__thread_local_fix double *A_rhs_array = NULL;
//up and down shade region
__thread_local_fix double *A_us_array;
__thread_local_fix double *A_ds_array;
//left and right shade  region
__thread_local_fix double *A_ls_array;
__thread_local_fix double *A_rs_array;
//front and back shade region
__thread_local_fix double *A_fs_array;
__thread_local_fix double *A_bs_array;
//Record the logical information of this time executed task
__thread_local_fix int dest,src,i,j,k,s_iter,s_starti,s_endi,s_startj,s_endj,s_startk,s_endk;
//Record the logical information of the last time executed task.
__thread_local_fix int last_iter,last_starti,last_startj,last_startk;
extern double A[2][SIZEX][SIZEY][SIZEZ];

#define alpha_d 0.0876
#define beta_d  0.0765
#define A_lhs(i,j,k) A_lhs_array[(i)*(s_endj-s_startj)*(s_endk-s_startk)+(j)*(s_endk-s_startk)+(k)]
#define A_rhs_unreuse(i,j,k) A_rhs_array[(i)*(s_endj-s_startj+2)*(s_endk-s_startk+2)+(j)*(s_endk-s_startk+2)+(k)]
#define A_rhs(i,j,k) A_rhs_array[(i)*(s_endj-s_startj)*(s_endk-s_startk)+(j)*(s_endk-s_startk)+(k)]
#define A_us(i,j) A_us_array[(i)*(s_endk-s_startk)+(j)]
#define A_ds(i,j) A_ds_array[(i)*(s_endk-s_startk)+(j)]
#define A_ls(i,j) A_ls_array[(i)*(s_endk-s_startk)+(j)]
#define A_rs(i,j) A_rs_array[(i)*(s_endk-s_startk)+(j)]
#define A_fs(i,j) A_fs_array[(i)*(s_endj-s_startj)+(j)]
#define A_bs(i,j) A_bs_array[(i)*(s_endj-s_startj)+(j)]

void stencil_core(targs* args)
{  
   //0. assumptions       
   //data region: [starti-1,endi] [startj-1,endj] [startk-1,endk]  
   //compute region: [starti, endi) [startj,endj) [startk,endk)
   int lhs_size;
   int rhs_size;
   
   s_iter = args->iter;
   dest = (s_iter) % 2 ;
   src = 1 - dest;
   s_starti = args->starti;
   s_endi = args->endi;
   s_startj = args->startj;
   s_endj = args->endj;
   s_startk = args->startk;
   s_endk = args->endk;
  
   lhs_size = (s_endi-s_starti)*(s_endj-s_startj)*(s_endk-s_startk);
   rhs_size = (s_endi-s_starti+2)*(s_endj-s_startj+2)*(s_endk-s_startk+2);
   if(A_lhs_array == NULL || A_rhs_array == NULL )
   {
      A_lhs_array = (double *)ldm_malloc(lhs_size * sizeof(double));
      A_rhs_array = (double *)ldm_malloc(rhs_size * sizeof(double));
   }
if((s_iter-1) == last_iter && s_starti == last_starti && s_startj == last_startj && s_startk == last_startk)
{
  //task hit, direct vertical reuse, non-boundary dara already in LDM!
  int stride,len,bsize;
  A_us_array = &A_rhs_array[lhs_size];
  A_ds_array = &A_rhs_array[lhs_size+(s_endj-s_startj)*(s_endk-s_startk)];
  A_ls_array = &A_rhs_array[lhs_size+(s_endj-s_startj)*(s_endk-s_startk)*2];
  A_rs_array = &A_rhs_array[lhs_size+(s_endj-s_startj)*(s_endk-s_startk)*2+(s_endi-s_starti)*(s_endk-s_startk)];
  A_fs_array = &A_rhs_array[lhs_size+(s_endj-s_startj)*(s_endk-s_startk)*2+(s_endi-s_starti)*(s_endk-s_startk)*2];
  A_bs_array = &A_rhs_array[lhs_size+(s_endj-s_startj)*(s_endk-s_startk)*2+(s_endi-s_starti)*(s_endk-s_startk)*2+(s_endi-s_starti)*(s_endj-s_startj)];

  get_reply = 0;
  //1.1 update up and down shade data
  stride = (SIZEZ-(s_endk-s_startk))*8;
  len = (s_endj-s_startj)*(s_endk-s_startk)*8;
  bsize = (s_endk-s_startk)*8;
  athread_get(PE_MODE,&A[src][s_starti-1][s_startj][s_startk]
                 , &A_us(0,0)
                 , len, &get_reply
                 , 0, stride, bsize);
  athread_get(PE_MODE,&A[src][s_endi][s_startj][s_startk]
                 , &A_ds(0,0)
                 , len, &get_reply
                 , 0, stride, bsize);
  //1.2 update left and right shade data
  stride = ((SIZEY-1)*SIZEZ+(SIZEZ-(s_endk-s_startk)))*8;
  len = (s_endi-s_starti)*(s_endk-s_startk)*8;
  bsize = (s_endk-s_startk)*8;
  athread_get(PE_MODE,&A[src][s_starti][s_startj-1][s_startk]
                   , &A_ls(0,0)
                   , len, &get_reply
                   , 0, stride, bsize);
  athread_get(PE_MODE,&A[src][s_starti][s_endj][s_startk]
                   , &A_rs(0,0)
                   , len,&get_reply
                   , 0, stride, bsize);    
  //2.compute core data
  for(i = 1; i < s_endi-s_starti-1; i++)
    for(j = 1; j < s_endj-s_startj-1; j++)
      for(k = 1; k < s_endk-s_startk-1; k++) 
         A_rhs(i,j,k) = alpha_d * (A_lhs(i,j,k)) + \
                     beta_d * (A_lhs(i-1,j,k) + A_lhs(i+1,j,k) +\
                     A_lhs(i,j-1,k) + A_lhs(i,j+1,k)+\
                     A_lhs(i,j,k-1) + A_lhs(i,j,k+1));
  //3.wait for up down left right data finished 
  while(get_reply != (2+2));
  
  get_reply = 0;
  //4. update front and back shade data
  stride = (SIZEZ-1)*8;
  len = (s_endj-s_startj)*8;
  bsize = 8;
  for(i = s_starti; i < s_endi; ++i )
  {
     athread_get(PE_MODE,&A[src][i][s_startj][s_startk-1]
                   , &A_fs(i-s_starti,0)
                   , len,&get_reply
                   , 0, stride, bsize);
     athread_get(PE_MODE,&A[src][i][s_startj][s_endk]
                   , &A_bs(i-s_starti,0)
                   , len,&get_reply
                   , 0, stride, bsize); 
  }
  //5.compute 4 line , 4 face.
  //5.1 compute 4 line:compute up and down line of k dimensional change
  for(k = 1; k < s_endk-s_startk-1; ++k)
  {
    //left face
    A_rhs(0,0,k) = alpha_d * (A_lhs(0,0,k)) + \
                     beta_d * (A_us(0,k) + A_lhs(1,0,k) +\
                     A_ls(0,k) + A_lhs(0,1,k)+\
                     A_lhs(0,0,k-1) + A_lhs(0,0,k+1));
    A_rhs(s_endi-s_starti-1,0,k) = alpha_d * (A_lhs(s_endi-s_starti-1,0,k)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,0,k) + A_ds(0,k) +\
                     A_ls(s_endi-s_starti-1,k) + A_lhs(s_endi-s_starti-1,1,k)+\
                     A_lhs(s_endi-s_starti-1,0,k-1)+ A_lhs(s_endi-s_starti-1,0,k+1));
   //right face 
    A_rhs(0,s_endj-s_startj-1,k) = alpha_d * (A_lhs(0,s_endj-s_startj-1,k)) + \
                     beta_d * (A_us(s_endj-s_startj-1,k) + A_lhs(1,s_endj-s_startj-1,k) +\
                     A_lhs(0,s_endj-s_startj-2,k) + A_rs(0,k)+\
                     A_lhs(0,s_endj-s_startj-1,k-1) + A_lhs(0,s_endj-s_startj-1,k+1));
    A_rhs(s_endi-s_starti-1,s_endj-s_startj-1,k) = alpha_d * (A_lhs(s_endi-s_starti-1,s_endj-s_startj-1,k)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,s_endj-s_startj-1,k) + A_ds(s_endj-s_startj-1,k) +\
                     A_lhs(s_endi-s_starti-1,s_endj-s_startj-2,k) + A_rs(s_endi-s_starti-1,k)+\
                     A_lhs(s_endi-s_starti-1,s_endj-s_startj-1,k-1) + A_lhs(s_endi-s_starti-1,s_endj-s_startj-1,k+1));
  }  
  //5.2 compute 2 face:compute left and right face
  for(i = 1; i < s_endi-s_starti-1; ++i)
    for(k = 1; k < s_endk-s_startk-1; ++k)
    {
       A_rhs(i,0,k) = alpha_d * (A_lhs(i,0,k)) + \
                     beta_d * (A_lhs(i-1,0,k) + A_lhs(i+1,0,k) +\
                     A_ls(i,k) + A_lhs(i,1,k)+\
                     A_lhs(i,0,k-1) + A_lhs(i,0,k+1));
       A_rhs(i,s_endj-s_startj-1,k) = alpha_d * (A_lhs(i,s_endj-s_startj-1,k)) + \
                     beta_d * (A_lhs(i-1,s_endj-s_startj-1,k) + A_lhs(i+1,s_endj-s_startj-1,k) +\
                     A_lhs(i,s_endj-s_startj-2,k) + A_rs(i,k)+\
                     A_lhs(i,s_endj-s_startj-1,k-1) + A_lhs(i,s_endj-s_startj-1,k+1));
    }
   //5.3compute 2 face: compute up and down face
   for(j = 1; j < s_endj-s_startj-1; ++j)
    for(k = 1; k < s_endk-s_startk-1; ++k)
    {
      A_rhs(0,j,k) = alpha_d * (A_lhs(0,j,k)) + \
                     beta_d * (A_us(j,k) + A_lhs(1,j,k) +\
                     A_lhs(0,j-1,k) + A_lhs(0,j+1,k)+\
                     A_lhs(0,j,k-1) + A_lhs(0,j,k+1));
      A_rhs(s_endi-s_starti-1,j,k) =  alpha_d * (A_lhs(s_endi-s_starti-1,j,k)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,j,k) + A_ds(j,k) +\
                     A_lhs(s_endi-s_starti-1,j-1,k) + A_lhs(s_endi-s_starti-1,j+1,k)+\
                     A_lhs(s_endi-s_starti-1,j,k-1) + A_lhs(s_endi-s_starti-1,j,k+1));
   }
  //6. wait for front back shade data  finished
  while(get_reply != ((s_endi-s_starti)*2));
  //7. compute 8 corner point
  //7.1 compute up_left_front_point
  A_rhs(0,0,0) = alpha_d * (A_lhs(0,0,0)) + \
                   beta_d * (A_us(0,0) + A_lhs(1,0,0) +\
                   A_ls(0,0) + A_lhs(0,1,0)+\
                   A_fs(0,0) + A_lhs(0,0,1)); 
  //7.2 compute up_left_back_point
  A_rhs(0,0,s_endk-s_startk-1) = alpha_d * (A_lhs(0,0,s_endk-s_startk-1)) + \
                   beta_d * (A_us(0,s_endk-s_startk-1) +\
                   A_lhs(1,0,s_endk-s_startk-1) +\
                   A_ls(0,s_endk-s_startk-1) + A_lhs(0,1,s_endk-s_startk-1)+\
                   A_lhs(0,0,s_endk-s_startk-2) + A_bs(0,0));
  
  //7.3 compute up_right_front_point
  A_rhs(0,s_endj-s_startj-1,0) = alpha_d * (A_lhs(0,s_endj-s_startj-1,0)) + \
                     beta_d * (A_us(s_endj-s_startj-1,0) +\
                     A_lhs(1,s_endj-s_startj-1,0) +\
                     A_lhs(0,s_endi-s_starti-2,0) + A_rs(0,0)+\
                     A_fs(0,s_endj-s_startj-1) + A_lhs(0,s_endj-s_startj-1,1));
  //7.4 compute up_right_back_point
  A_rhs(0,s_endj-s_startj-1,s_endk-s_startk-1) = alpha_d * \
                    (A_lhs(0,s_endj-s_startj-1,s_endk-s_startk-1)) + \
                     beta_d * (A_us(s_endj-s_startj-1,s_endk-s_startk-1) + \
                     A_lhs(1,s_endj-s_startj-1,s_endk-s_startk-1) + \
                     A_lhs(0,s_endj-s_startj-2,s_endk-s_startk-1) + \
                     A_rs(0,s_endk-s_startk-1)+\
                     A_lhs(0,s_endj-s_startj-1,s_endk-s_startk-2) + \
                     A_bs(0,s_endj-s_startj-1));  
  //7.5 compute down_left_front_point
   A_rhs(s_endi-s_starti-1,0,0) =   alpha_d * (A_lhs(s_endi-s_starti-1,0,0)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,0,0) + A_ds(0,0) +\
                     A_ls(s_endi-s_starti-1,0) + A_lhs(s_endi-s_starti-1,1,0)+\
                     A_fs(s_endi-s_starti-1,0) + A_lhs(s_endi-s_starti-1,0,1));  
  //7.6 compute down_left_back_point
   A_rhs(s_endi-s_starti-1,0,s_endk-s_startk-1) = alpha_d * \
                     (A_lhs(s_endi-s_starti-1,0,s_endk-s_startk-1)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,0,s_endk-s_startk-1) + \
                     A_ds(0,s_endk-s_startk-1) +\
                     A_ls(s_endi-s_starti-1,s_endk-s_startk-1) +\
                     A_lhs(s_endi-s_starti-1,1,s_endk-s_startk-1)+\
                     A_lhs(s_endi-s_starti-1,0,s_endk-s_startk-2) +\
                     A_bs(s_endi-s_starti-1,0));    
  //7.7 compute down_right_front_point
  A_rhs(s_endi-s_starti-1,s_endj-s_startj-1,0) = alpha_d * \
                     (A_lhs(s_endi-s_starti-1,s_endj-s_startj-1,0)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,s_endj-s_startj-1,0) + \
                     A_ds(s_endj-s_startj-1,0) + \
                     A_lhs(s_endi-s_starti-1,s_endj-s_startj-2,0) + \
                     A_rs(s_endi-s_starti-1,0)+ \
                     A_fs(s_endi-s_starti-1,s_endj-s_startj-1) + \
                     A_lhs(s_endi-s_starti-1,s_endj-s_startj-1,1)); 
  //7.8  compute down_right_back_point
  A_rhs(s_endi-s_starti-1,s_endj-s_startj-1,s_endk-s_startk-1) = alpha_d * \
                     (A_lhs(s_endi-s_starti-1,s_endj-s_startj-1,s_endk-s_startk-1)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,s_endj-s_startj-1,s_endk-s_startk-1) + \
                     A_ds(s_endj-s_startj-1,s_endk-s_startk-1) + \
                     A_lhs(s_endi-s_starti-1,s_endj-s_startj-2,s_endk-s_startk-1) + \
                     A_rs(s_endi-s_starti-1,s_endk-s_startk-1)+ \
                     A_lhs(s_endi-s_starti-1,s_endj-s_startj-1,s_endk-s_startk-2) + \
                     A_bs(s_endi-s_starti-1,s_endj-s_startj-1)); 
 //8 compute 8 line data
  //8.1 compute up and down line of j dimensional  change 
  for(j = 1;j < s_endj-s_startj-1; ++j)
  {
  //front face
    A_rhs(0,j,0) = alpha_d * (A_lhs(0,j,0)) + \
                     beta_d * (A_us(j,0) + A_lhs(1,j,0) +\
                     A_lhs(0,j-1,0) + A_lhs(0,j+1,0)+\
                     A_fs(0,j) + A_lhs(0,j,1));
    A_rhs(s_endi-s_starti-1,j,0) = alpha_d * (A_lhs(s_endi-s_starti-1,j,0)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,j,0) + A_ds(j,0) +\
                     A_lhs(s_endi-s_starti-1,j-1,0) + A_lhs(s_endi-s_starti-1,j+1,0)+\
                     A_fs(s_endi-s_starti-1,j) + A_lhs(s_endi-s_starti-1,j,1));
  //back face
    A_rhs(0,j,s_endk-s_startk-1) = alpha_d * (A_lhs(0,j,s_endk-s_startk-1)) + \
                     beta_d * (A_us(j,s_endk-s_startk-1) + A_lhs(1,j,s_endk-s_startk-1) +\
                     A_lhs(0,j-1,s_endk-s_startk-1) + A_lhs(0,j+1,s_endk-s_startk-1)+\
                     A_lhs(0,j,s_endk-s_startk-2) + A_bs(0,j));
    A_rhs(s_endi-s_starti-1,j,s_endk-s_startk-1) =  alpha_d * (A_lhs(s_endi-s_starti-1,j,s_endk-s_startk-1)) + \
                     beta_d * (A_lhs(s_endi-s_starti-2,j,s_endk-s_startk-1) + A_ds(j,s_endk-s_startk-1) +\
                     A_lhs(s_endi-s_starti-1,j-1,s_endk-s_startk-1) + A_lhs(s_endi-s_starti-1,j+1,s_endk-s_startk-1)+\
                     A_lhs(s_endi-s_starti-1,j,s_endk-s_startk-2) + A_bs(s_endi-s_starti-1,j));
  }
  //8.2 compute left and right line of i dimensional change
  for(i = 1;i < s_endi-s_starti-1; ++i)
  {
   //front face
    A_rhs(i,0,0) = alpha_d * (A_lhs(i,0,0)) + \
                     beta_d * (A_lhs(i-1,0,0) + A_lhs(i+1,0,0) +\
                     A_ls(i,0) + A_lhs(i,1,0)+\
                     A_fs(i,0) + A_lhs(i,0,1));
    A_rhs(i,s_endj-s_startj-1,0) = alpha_d * (A_lhs(i,s_endj-s_startj-1,0)) + \
                     beta_d * (A_lhs(i-1,s_endj-s_startj-1,0) + A_lhs(i+1,s_endj-s_startj-1,0) +\
                     A_lhs(i,s_endj-s_startj-2,0) + A_rs(i,0)+\
                     A_fs(i,s_endj-s_startj-1) + A_lhs(i,s_endj-s_startj-1,1)); 
  //back face 
    A_rhs(i,0,s_endk-s_startk-1) = alpha_d * A_lhs(i,0,s_endk-s_startk-1) + \
                     beta_d * (A_lhs(i-1,0,s_endk-s_startk-1) + A_lhs(i+1,0,s_endk-s_startk-1) +\
                     A_ls(i,s_endk-s_startk-1) + A_lhs(i,1,s_endk-s_startk-1)+\
                     A_lhs(i,0,s_endk-s_startk-2) + A_bs(i,0)); 
    A_rhs(i,s_endj-s_startj-1,s_endk-s_startk-1) = alpha_d * (A_lhs(i,s_endj-s_startj-1,s_endk-s_startk-1)) + \
                     beta_d * (A_lhs(i-1,s_endj-s_startj-1,s_endk-s_startk-1) + A_lhs(i+1,s_endj-s_startj-1,s_endk-s_startk-1) +\
                     A_lhs(i,s_endj-s_startj-2,s_endk-s_startk-1) + A_rs(i,s_endk-s_startk-1)+\
                     A_lhs(i,s_endj-s_startj-1,s_endk-s_startk-2) + A_bs(i,s_endj-s_startj-1));
  }
  //9 compute 2 face
  //9.1 compute front and back face 
  for(i = 1; i < s_endi-s_starti-1; ++i)
    for(j = 1; j < s_endj-s_startj-1; ++j)
    {
       A_rhs(i,j,0) = alpha_d * (A_lhs(i,j,0)) + \
                     beta_d * (A_lhs(i-1,j,0) + A_lhs(i+1,j,0) +\
                     A_lhs(i,j-1,0) + A_lhs(i,j+1,0)+\
                     A_fs(i,j) + A_lhs(i,j,1)); 
       A_rhs(i,j,s_endk-s_startk-1) = alpha_d * (A_lhs(i,j,s_endk-s_startk-1)) + \
                     beta_d * (A_lhs(i-1,j,s_endk-s_startk-1) + A_lhs(i+1,j,s_endk-s_startk-1) +\
                     A_lhs(i,j-1,s_endk-s_startk-1) + A_lhs(i,j+1,s_endk-s_startk-1)+\
                     A_lhs(i,j,s_endk-s_startk-2) + A_bs(i,j));  
    }
//10 write back
   put_reply=0;
   for(i =  s_starti; i < s_endi; ++i)
   {   
     int len,bsize,stride;
     stride = (SIZEZ-(s_endk-s_startk))*8;
     bsize = (s_endk-s_startk)*8;
     len = bsize*(s_endj-s_startj);
     athread_put(PE_MODE,&A_rhs(i-s_starti,0,0)
                     , &A[dest][i][s_startj][s_startk]
                     ,len,&put_reply,stride,bsize);
   }   
   while(put_reply!=(s_endi-s_starti)); 
}else{
// task miss, no direct vertical reuse, all ref data from main memory
//1 update core data
  get_reply=0; 
  for(i = s_starti-1; i < s_endi+1; i++)
  {
       int len, bsize, stride;
       bsize=(s_endk-s_startk+2)*8;
       stride=(SIZEZ-(s_endk-s_startk+2))*8;  //stride: from blocktail to blockhead 
       len= bsize*(s_endj-s_startj+2);
       athread_get(PE_MODE,&A[src][i][s_startj-1][s_startk-1]
                 , &A_rhs_unreuse(i-s_starti+1,0,0)
                 , len, &get_reply
                 , 0, stride, bsize);
  }
  //2.wait for core data finished    
  while(get_reply!=(s_endi-s_starti+2));
  //3.compute core data
  for(i = 1; i < s_endi-s_starti+1; i++)
     for(j = 1; j < s_endj-s_startj+1; j++)
        for(k = 1; k < s_endk-s_startk+1; k++) 
            A_lhs(i-1,j-1,k-1) = alpha_d * (A_rhs_unreuse(i,j,k)) + \
                     beta_d * (A_rhs_unreuse((i-1),j,k) + A_rhs_unreuse((i+1),j,k) +\
                     A_rhs_unreuse(i,(j-1),k) + A_rhs_unreuse(i,(j+1),k)+\
                     A_rhs_unreuse(i,j,(k-1)) + A_rhs_unreuse(i,j,(k+1)));     
  //4.write back
   put_reply=0;
   for(i =  s_starti; i < s_endi; ++i)
   {   
     int len,bsize,stride;
     stride = (SIZEZ-(s_endk-s_startk))*8;
     bsize = (s_endk-s_startk)*8;
     len = bsize*(s_endj-s_startj);
     athread_put(PE_MODE,&A_lhs(i-s_starti,0,0)
                     , &A[dest][i][s_startj][s_startk]
                     ,len,&put_reply,stride,bsize);
   }   
   while(put_reply!=(s_endi-s_starti)); 
   last_iter = s_iter;
   last_starti = s_starti;
   last_startj = s_startj;
   last_startk = s_startk;
}
//5.free A_slave_array AND return to scheduler loop.
   ldm_free(A_lhs_array,lhs_size*sizeof(double));
   ldm_free(A_rhs_array,rhs_size*sizeof(double));
   
}
